Kramers escape driven by fractional Brownian motion 
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We investigate the Kramers escape from a potential well of a test particle driven by fractional 
Gaussian noise with Hurst exponent < -ff < 1. From a numerical analysis we demonstrate the 
exponential distribution of escape times from the well and analyze in detail the dependence of the 
mean escape time as function of H and the particle diffusivity D. We observe different behavior 
for the subdiffusive (antipersistent) and superdiffusive (persistent) domains. In particular we find 
that the escape becomes increasingly faster for decreasing values of H, consistent with previous 
findings on the first passage behavior. Approximate analytical calculations are shown to support 
the numerically observed dependencies. 

PACS numbers: 05.40.Fb,02.50.Ey 



I. INTRODUCTION 

Anomalous diffusion is characterized by a deviation 
from the classical linear time dependence {x^it)) ~ t of 
the mean squared displacement. Such anomalies range 
from ultraslow transport {x^) ~ loS^^ discovered in 
Sinai diffusion or in iterated maps i up to cubic diffu- 
sion {x{t)) ~ in random walk processes with correlated 
jump lengths [3] or the relative coordinate of two particles 
encountered in turbulent Richardson flow [4, 5] . Here we 
are interested in anomalous diffusion of the power-law 
type 0,0 



''{t)) = 2Dt 



2H 



(1) 



where H is the Hurst exponent and D the generalized 
diffusion coefficient of dimension [D] = cm^/sec^^. De- 
pending on the magnitude of H we observe subdiffu- 
sion (0 < iJ < 1/2) or superdiffusion (1/2 < H < 1). 
The limits H — 1/2 and H = 1 correspond to ordinary 
Brownian diffusion or ballistic motion, respectively. For 
one-particle motion ballistic transport is the upper limit 
of spreading when the particle has a finite maximum ve- 
locity. 

Anomalous diffusion of the power law form ([T]) is ob- 
served in a multitude of systems. In particular, sub- 
diffusion was found for the motion of charge carriers in 
amorphous semiconductors , the spreading of tracer 
molecules in subsurface hydrology [lo| . diffusion on ran- 
dom site percolation clusters [ll| as well as the motion 
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of tracers in the crowded environment of biological cells 



[12j or in reconstituted biological systems [13|, among 
many others. Examples for superdiffusion include active 
motion in biological cells [wj , tracer spreading in layered 
velocity fields |15| . turbulent rotating flows [16], or in 
bulk mediated surface exchange jl7| . 

Apart from numerical approaches there exist two 
prominent analytical models for such anomalous diffu- 
sion: One is the continuous time random walk (CTRW) 
model [1, [Tsj in which each jump is characterized by a 
variable jump length and waiting time drawn from asso- 
ciated probability densities. CTRW theory includes (i) 
subdiffusion when the variance of jump lengths is finite 
but the waiting times have an infinite characteristic time; 
(ii) Levy flights when the mean waiting time is finite but 
the jump length variance diverges; and (iii) Levy walks in 
which waiting times and jump lengths are coupled, pro- 
ducing sub-ballistic superdiffusion with finite variance. 
The escape over a potential barrier for subdiffusion and 
Levy flights was studied recently [l9l - l22j [. 

The second model is fractional Brownian motion 
(FBM). It was originally described by Kolmogorov [23| 
and reintroduced by Mandelbrot and van Ness I24| . FBM 
is a self-similar Gaussian process with stationary incre- 
ments [1^, [2^. The FBM mean squared displacement 
follows Eq. ((TJ, and the Hurst exponent H of the frac- 
tional Gaussian noise varies in the full range < H < 1. 
Uncorrelated, regular Brownian motion corresponds to 
H = 1/2. For < H < 1/2 the prefactor in the noise 
autocorrelation is negative, rendering the associated an- 
tipersistent process subdiffusive. That means that a step 
in one direction is likely followed by a step in the other 
direction. Conversely, in the case 1/2 < H < 1 the mo- 
tion is persistent, effecting sub-ballistic superdiffusion in 
which successive steps tend to point in the same direc- 
tion. FBM is used to model a variety of processes includ- 
ing monomer diffusion in a polymer chain [27| , single flle 
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Figure 1: Harmonic potential V{x) — with cutoff at 

X = y/2 used in the numerical analysis of the FBM Kramers 
escape. The potential barrier height is AV^ = 1. Dimension- 
less units. 



diffusion [28'|, diffusion of biopolymers in the crowded en- 
vironment inside biological cells ^9], long term storage 
capacity in reservoirs [s^], climate fluctuations [3l|, eco- 
nophysics [32|, and teletraffic [ssj . 

Despite its wide use FBM is not completely under- 
stood. Thus the general incorporation of non-trivial 
boundary conditions is unattained, in particular, the 
first passage behavior is solved analytically solely on a 
semi- infinite domain [sj]- Notably the method of im- 
ages does not apply to solve boundary value problems 
for FBM. Similarly the associated fractional Langevin 
equation driven by fractional Gaussian noise was recently 
discovered to exhibit critical dynamical behavior [35]. 

Here we study the generalization for FBM of the 
Kramers escape from a potential well across a finite bar- 
rier, as illustrated in Fig. [T] This problem is relevant, 
for instance, for single file diffusion in external poten- 
tials [IBl , the dissociation dynamics of biopolymers from 
a bound state in FBM models for particle diffusion un- 
der molecular crowding conditions [29| or bulk chemical 
reactions of larger particles under superdense conditions. 
We note that a similar problem was treated for correlated 
Gaussian noise and for fractional Langevin equation 
motion in the case when the fluctuation dissipation the- 
orem applies dl]. We here study the important case of 
external fluctuations, that is, for systems which do not 
obey the fluctuation dissipation theorem [39| . 

In the regular Kramers theory [40l - |42j for the escape of 
a Brownian particle across a potential barrier in the high 
barrier limit AV ^ ksO, where ksd denotes thermal 
energy, the probability density of the first escape from 
the well follows an exponential decay, 



p{t) ~ exp ( - ^ 



(2) 



eigenvalue |4CH42| | . In Eq. ([2]) the characteristic escape 
time T is proportional to the Arrhenius factor of the bar- 
rier height AV, 



(3) 



In what follows we demonstrate from simulations and 
analytical considerations that the exponential decay ([2]) 
is preserved in FBM processes due to the stationary 
nature of FBM, while the activation pattern ^ becomes 
explicitly dependent on the Hurst exponent. This H- 
dependence is different for the antipersistent and persist- 
ent cases. Remarkably slow diffusion leads to fast escape, 
that is, the lower the value of H is chosen the faster the 
escape from the potential well becomes. This observa- 
tion is consistent with the first passage behavior of FBM 
that is known analytically, and analysed numerically in 
the Appendix. 

We first investigate FBM driven Kramers escape by 
numerical integration of the Langevin equation subject 
to fractional Gaussian noise in Sec. HIl In particular, we 
analyze the distribution of escape times and the depend- 
ence of the mean escape time on the Hurst exponent H 
and the noise strength D. In Sec. IIIII we develop an ap- 
proximate analytical approach to the barrier crossing for 
FBM, before drawing our conclusions in Sec. IIVI In the 
Appendices we describe the numerical algorithms used to 
generate antipersistent and persistent FBM, and we val- 
idate in detail that these truthfully produce FBM. We 
also briefly discuss the consistency of our results for the 
case of a potential well, that is finite on both sides. 



II. NUMERICAL ANALYSIS 

In this Section we set up the Langevin description of 
FBM for external Gaussian noise and present extensive 
simulations results for the barrier crossing behavior. 



A. Langevin equation with fractional Gaussian 
noise 

We employ the overdamped Langevin equation for the 
position variable x{t) in the presence of an external po- 
tential V{x), 



dx{t) 
dt 



1 dV{x) 
dx 



7717 



(4) 



where m is the particle mass, 7 the friction constant, 
^H{t) is the fractional Gaussian noise, and D is its in- 
tensity. The chosen initial condition is a;(0) = 0. To 
study the activated escape from a potential well, in what 
follows we use an harmonic potential of the form 



This corresponds to the relaxation mode of the lowest 



V{x)^ 



—00 < X < \f2 



-00, X > \/2 



(5) 



3 



with a truncation at positive x = \f2, compare Fig. [T] 
We note that we compared our simulations for the po- 
tential ([5]) to the escape from an harmonic potential with 
symmetric truncation, 

(cx), — oo < X < -~\pl 
§a;2, <x<V2 , (6) 

—oo, X > \/2 

finding qualitative agreement with the results reported 
herein with respect to the dependence of the distribution 
of escape times and the dependence of the mean escape 
time on Hurst exponent and noise strength, see the Ap- 
pendix. 

In continuous time the fractional Gaussian noise (t) 
is understood as a derivative of the FBM [H [2^. This 
is a stationary Gaussian process with an autocorrelation 
function that in the long time limit decays as 

{^HmH{t))-^2H{2H-l)t^"-^ (7) 

ioT < H < 1, H 1/2. Note that in the antipersistent 
case, < H < 1/2, the autocorrelation function of the 
fractional Gaussian noise is negative at long times. At 
H ^ 1/2 we have a delta-correlated white noise. In a dis- 
crete time approximation used in numerical simulations 
below the autocorrelation function of the noise reads (2^ 

{^HmHin)) = [{n + If" - 2n'" + {n - if"] . (8) 

The continuum approximation ([T]) is obtained from 
Eq. ([SJ in the limit of large n and identifying n ^ t. 
In what follows in analytical calculations and numerical 
simulations we use the PDF of the fractional Gaussian 
noise 




with variance 2. 

Replacing x — (7717/0)^^^2; and t — > {m'f/a)t we 
pass to reduced variables: 

^^-X + V^^Hit). (10) 

The time-discretized version of Eq. ([TU)) acquires the form 

Xn+l - Xn = -XnSt + VOSt^Hin), (11) 

where St is a finite time step. 

We applied the methods described in Refs. [i^ and [i^l 
for simulating fractional Gaussian noise with H < 1/2 
and H > 1/2, respectively, as detailed in Appendix 1X1 In 
the simulations the Hurst index H was varied within the 
range [0.1,0.85], whereas the noise intensity D covered 
values from 1/6 to 1/2. Correspondingly, the escape 
time was varying in a range covering three orders of mag- 
nitude. 
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Figure 2: Top panel: Probability density function (PDF) of 
the escape time nicely demonstrating the exponential decay 
given by Eq. ([2]). The main plot depicts the persistent case 
{H = 0.55, 0.65, 0.75), while in the inset we show the anti- 
persistent case {H — 0.1, 0.2, and 0.3). Here we used the 
following simulations parameters: in the antipersistent case 
the time increment is St — 0.001, the number of samples 
A'^stat ~ 10'', the number of data points per sample A^max 
varied from 2" 8.2 x 10^ to 2^^ 2.1 x 10®, and finally 
the noise strength D = 0.25; in the persistent case we used 
St = 0.001, iVstat = 20000, A^max = 2^° « 10'', and D = 0.25. 
Bottom panel: PDF of the escape time in a log-log repres- 
entation. The decay curves have a common envelope 1 / (et) 
depicted by the straight grey line, see text. 



B. Numerical results for FBM Kramers escape 

In our simulations we follow the motion of the test 
particle governed by the discrete Langevin equation (jlip 
in the harmonic potential with one-sided truncation, 
Eq. (O. Once the particle crosses the point x = \/2 
it is removed, and the next particle started. This setup 
is depicted in Fig. [1] 

We first focus on the probability density function 
(PDF) of the first escape time from the potential well. 
In Fig. [5] we demonstrate that, in analogy to the classical 
case (_ff = 1/2) the probability density function (PDF) 
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Figure 3: Top panel: Mean escape time T as function of 
inverse noise intensity 1/D in logarithmic scale. The simula- 
tions results are shown for Hurst exponents H = 0.85, 0.75, 
0.65, 0.55, 0.5, 0.4, 0.3, 0.2, and 0.1 (top to bottom). The solid 
lines represent a linear fit. Bottom panel: Fitting coefficients 
a{H) and b{H) from Eq. ([12}. The inset shows b{H) for the 
antipersistent case at higher resolution. The symbols repres- 
ent the values of a and b from the simulations, while the solid 
lines show the respective fits with Eqs. H13p to (|16p . Simula- 
tions parameters: in the antipersistent case we used the time 
increment 5t = 0.001, number of samples A^stat ~ 10^, and 
the data points per sample A'max varied from 2^^ ~ 8.2 x 10"^ 
to 2^^ ~ 2.1 X lO''; in the persistent case 5t ranges from 0.01 
to 0.001, Nst!,t = 10^, and iVmax varied from 2" 8.2 x 10^ 
to 2^** « 2.6 X 10^ 



of the first escape time decays exponentially with time, 
see Eq. This exponential decay is observed nicely in 
the simulations data over the entire range of the Hurst 
exponent. In the double-logarithmic plot in the bottom 
panel of Fig. [2] one can see a common envelope of the 
curves for all values of H. Indeed the shoulders of the in- 
dividual exponential PDFs are located at points in time 
where t = T, i.e., where the value of the PDFs is ex- 
actly l/(eT). This is the straight line plotted in Fig. [21 
showing good agreement, with a slight underestimation 
for persistent Hurst exponents. 



1 2 3 4 5 6 

1/D 

Figure 4: Mean escape time T as function of inverse noise 
intensity 1/D in linear scale. Data points are the same as 
in Fig. [3l Again, the dots represent the simulations data, 
and the solid lines depict the exponential fit with the global 
fitting coefficients according to Eqs. H13I) to p6|) . From top 
to bottom: Hurst exponent H = 0.85 (V), H = 0.75 (x), 
H = 0.65 (□), H = 0.55 (A), H = 0.50 {+), H = 0.40 (V), 
H = 0.30 ((□), H = 0.20 (A), H = 0.10 (x). Inset: Too high 
noise intensity violating the high barrier assumption leads to 
deviations from the exponential behavior. Hurst exponent 
H = 0.5 (□) and H = 0.8 (x). 

In Fig. [3] we demonstrate that the mean escape time T 
follows an exponential behavior as function of the inverse 
noise intensity, in analogy to the classical Kramers 

case. We observe that in both persistent and antipersist- 
ent cases this functional dependence may be approxim- 
ated by a linear fit of the form 

lnT(D;H)=a{H) + ^-^, (12) 

where both fitting coefficients a and b are functions of 
the Hurst exponent H. These, in turn, show different 
behavior for antipersitence and persistence of the motion: 

(i) In the persistent case 1/2 < iJ < 1 both coefficients 
are linear functions of the Hurst exponent. We found 
empirically from best fits that 

a{H> 1/2) = ai+a2H, (13) 
b{H > 1/2) = bi+b2H, (14) 

where oi = -1.680, aa = 4.869, bi = 1.051, and 62 = 
—0.399. The good quality of this linear description is 
seen in the bottom panel of Fig. |31 

(ii) Contrasting this behavior, in the antipersistent 
case the coefficient a (H) is still well described by a lin- 
ear i?-dependence, while b {H) is well represented by a 
parabolic dependence: 

a{H <l/2) = di+a2H, (15) 
b{H< 1/2) = bi+b2H + b3H^. (16) 
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Figure 5: Mean escape time as a function of the Hurst ex- 
ponent for different noise intensities D. The sohd hnes cor- 
respond to the fits used in Fig. |3l converted according to 
Eqs. (fT7)l and (fTS)) . 



2 



^ 0.5 



- 
















- 


antipersistent 


— 










persistent 



0.05 0.1 0.15 0.2 0.25 



|H-1/2| 

Figure 6: Relative escape rate k (H) /k (H — 1/2) as a func- 
tion of \H — 1/2| for the persistent and antipersistent cases 
for D = 0.25. The solid parabolic and straight lines are ob- 
tained from fits to the numerical data according to Eqs. (|17|) 
and (|18p . respectively. A fine coincidence is observed. 



The best fit parameters are determined as oi = —3.019, 
02 = 7.296, 6i = 0.705, h = 1.490 and 63 = -2.281. 
Again, Fig. [3] demonstrates good agreement with this 
chosen _ff-dependence. In Fig. 2] we show the quality 
of these fits (solid curves) on a linear scale. Note the 
deviations from the exponential behavior when the noise 
intensity becomes too large [in our simulation for values 
D > 1]. In that case the high barrier limit is violated 
and the results obtained herein are no more applicable, 
in correspondence to regular Brownian barrier crossing 
behavior. 

The general agreement with the law (IT^ is excellent, 
keeping in mind that the error of the simulations data 
is of the magnitude of the points. Remarkably the char- 
acteristic escape time increases from low to high Hurst 
exponent. In other words, the less persistent motion 
shows the faster escape. This observation is consistent 
throughout our simulations. In particular this behavior 
is not qualitatively changed for a parabolic potential of 
the type © with symmetric cutoff. 

In Fig. [5] the mean escape time is reanalyzed as a func- 
tion of the Hurst exponent. In accordance with the res- 
ults presented in Fig. |3l there is a parabolic dependence 
of InT versus H in the antipersistent case (0 < H < 1/2), 

lnT = ci + C2H + d3H^, (17) 

where ci ^ ai + bi/ D, C2 ~ a-i + ^2/0, and C3 — b^/ D. 
In the persistent case 1/2 < < 1 the relation is linear, 
corresponding to 

\nT = ci+C2H, (18) 

where ci — ai+bi/D and C2 = 02+62/!^. The agreement 
with the fit function is favorable, and the continuation 
between antipersistent and persistent cases appears relat- 
ively smooth. The latter supports the good convergence 



of the simulations algorithms used in the antipersistent 
and persistent regimes (see Appendix A). At the same 
time the difference between the behaviors in the two re- 
gimes (persistent versus antipersistent) is quite distinct. 

Fig. ini shows an alternative way to represent the be- 
havior from Figs. [3] and [SI namely, in terms of the ratio 
k (H) /k {H = 1/2) of the escape rates (that is, the in- 
verse mean escape times) as function of the deviation 
\H — 1/2 1 from normal difi'usion at H — 1/2. The rates 
increase with decreasing Hurst exponent, i.e., the less 
persistent the motion is the higher becomes the corres- 
ponding rate. One can also see the difference between 
the parabolic dependence in the antipersistent case and 
the linear relation for persistent motion. 

Finally in Fig. [7] we explore the distribution of the res- 
ults for the mean escape time between different samples 
of only 60 trajectories. Again we see the increased escape 
time at higher Hurst exponent. We also clearly observe 
that the variation around the average values increases 
significantly for higher Hurst exponent. In particular the 
noise for the plotted case H = 0.3 is consistently smaller 
than for the Brownian limit H = 1/2. 



III. ANALYTICAL APPROACH TO FBM 
DRIVEN KRAMERS ESCAPE 

In this Section we derive analytical results for the es- 
cape behavior driven by fractional Gaussian noise. In 
particular we concentrate on the mean escape time and 
the autocorrelation function for FBM in an harmonic po- 
tential. We compare the results to the numerical findings 
from the preceding Section. 
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Figure 7: (Color online) Results for the mean escape time. 
Main graph: for three different values of the Hurst expo- 
nents {H = 0.3, 0.5, and 0.7) we show the average value and 
the standard deviation around that value for different inverse 
noise intensities 1/D. Note that Tesc is evaluated at 1/_D = 2, 
2.5, 3, . . .. at each H value. In the figure the results for dif- 
ferent H at each given value of 1/D are slightly shifted with 
respect to each other for illustrative purposes. Inset: values of 
the escape times for each individual out of 60 trajectories for 
H — 0.7. The grey line connects the averages. In the simula- 
tions we used the time increment St — 0.001 and the number 
of points per trajectory A'max varied from 2^^ ~ 2 x 10'^ to 
2^0 ^ 10«. 



A. Wilemski-Fixmann approximation 

The investigation of first passage times for non- 
Markovian processes has a long historyin mathemat- 
ical literature, for instance, see Refs. jS^, |45| - |47| . and 
appears in different fields of science, including chemical 
physics [isj . polymer physics [i^, [s^ and neuroscience 
[5l| . However, no general theory exists for such pro- 
cesses, and different approximations are used depending 
on whether the process is Gaussian or not, whether its 
trajectories are differentiable or not, etc. Our analyt- 
ical approach to the escape problem considered herein 
is based on a special case of the Wilemski-Fixmann ap- 
proximation (WFA) [4§] used in polymer physics [s^l • As 
shown in Ref. 52] the application of the WFA to a first 
passage problem corresponds to a renewal approxima- 
tion [53, [SJ] in which, however, the correct Green's func- 
tions of the original non-renewal processes are used. The 
WFA is essentially a first approximation in the perturb- 
ative series derived by Likthman and Marques [53 , while 
higher approximations lead to quite involved expressions. 

Our theoretical approach starts from the relation 

Gix,t\xo,0) = S{x-xo)S{t) 

+ [ F{x,t' ,xa,0)G{x,t\x,t')dt' , (19) 



ted at xq Sit time t = Q. Moreover F{x, t, xq, 0) represents 
the first passage time PDF to cross the distance |a; — a;o| 
during the time interval t, and G {x, t\x, t') is the condi- 
tional probability to be at x at time t, provided x was 
visited earlier at time t' . If the inequality xq ^ x holds 
the (5-terni can be omitted. For a continuous Markovian 
process Eq. is exact. Its meaning is that a particle, 
having started at xq at time and being at a site x at 
time might have visited x at some time t' before, depar- 
ted from X, and returned (53l.[5^. For the non-Mar kovian 
case Eq. neglects the correlations in the motion of 
the particle before and after the first passage through 
the point x. Such correlations lead to the dependence of 
the return probability (expressed through G{x,t\x,t'y) 
on the pre- history [52], and can be taken into account 
systematically in higher order approximations involving 
multi-point distribution functions (55| . The approxima- 
tion given by Eq. (|19p may become incorrect in the case 
of strongly correlated (persistent) processes. In that case 
our numerical results still show exponential first passage 
time behavior corresponding to a finite mean first pas- 
sage time, while the WFA breaks down, as will be shown 
below. 

To proceed recall that according to Bayes' formula, 
G{x,t\xQ,W) = P(x,i;a;o,0)/P(xo,0) and G{x,t\x,t') = 
P{x,t;x,t')/P{x,t'). Here P{x,t;x,0) and Pix,t) are 
the corresponding two- and one-point probability densit- 
ies. Eq. (|19p can therefore be rewritten in the form 

P{x,t; xo,0) = P{xo,0) 

X f F[x,t',x,,o f^lfl]p dt'. (20) 

Integration with respect to xq in Eq. ([20)1 leads to the 
expression 

P{x,t)^ f F{x,t')^^0^dt', (21) 



P(x,t') 



where 



/oo 
P[xa,Q)F[x,t\xo,t)dxo. 
-00 



(22) 



Thus, the first escape PDF F is obtained as an average 
over the initial distribution. 

In what follows we make use of the fact that in our 
numerical simulations the typical relaxation times for a 
particle in an harmonic potential well are much shorter 
than the typical mean escape times. Therefore the ran- 
dom process x{t) can be considered as stationary, that is, 
Pix,t) = Pst{x) a.nAP{x,t]x,t') = P{x,x,t-t'). Trans- 
ferring P from the left hand side to the right of Eq. (pij) 
we find 

F{x,t') ^^''L':'^7^'^ dt' 



1 



(23) 



This relation converts to an algebraic equation after 
Laplace transformation, 



where G{x,t\xo,0) is the conditional probability to find 
the particle at position x at time t, provided that it star- 



1 £,/ ~,P{x,x,s) 
- = F{x,s) 



PU^) 



(24) 
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Here we express the Laplace transform of a function f{t) 
as /(s) = f{t)exp{—st)dt. Since P{x,x,t — oo) — > 
Pg^j(x), we see that P{x,x,s — > 0) — > P^f.{x)/s, and for 
small s we may expand P{x, x, s) in the form 



P{x, X, s) 



+ A{x) + Ois), 



(25) 



where we use the abbreviation 



Aix) 



lim 



P{x, X, s) 



P{x,x,t) - Pl{x) 



dt. 



(26) 



After inserting Eq. (1^ into Eq. ([M)) we get 

Pl{x) ^ PUx) 



F{x,s) 



sPix,x,s) P!tix)+A{x)s 

p^d^y 



(27) 



Thus, with the use of Eq. ([25]), we find 



T{x) 



as 



A{x) 



Pli^) 



P{x,x,t) 

PU^) 



1 



dt. (28) 



We will use this result below. 

Before proceeding two remarks are in order: First, we 
note that in the theory developed here we use the en- 
semble average over initial values xo, while in the simu- 
lations we use xo = for all trajectories. Nevertheless, 
we can employ Eq. (|27p since typically the relaxation 
time is much shorter than the mean escape time and, 
therefore, the system quickly converges to the stationary 
state, which is independent of the initial condition. And 
second, when writing Eq. (|25p we implicitly assume that 
the mean escape time exists. This is in accordance with 
the numerical observation that the escape time PDF has 
the simple exponential form 



B. Mean escape time for Gaussian processes 

To proceed we exploit the Gaussian property of FBM 
processes. We recall the expressions for one- and two- 
point Gaussian PDFs, namely, 



Pst{x) 



1 



V2t 



: exp 



X 

2a^ 



(29) 



where = (a;^)^^ is the variance in the stationary state 
of a particle in an harmonic potential well. Moreover 



Pix,y,t) = 



X exp 



2(1 -.g2) 



2gxy 



(30) 



where g(£) is the normalized autocorrelation function in 
the stationary state. 



{x{t)v{t + T))^ 



Thus, within our approximation 
1 



P{x, x, t) 



2^2 v/r^ 

and we obtain the mean time 



exp 



a2(l + .g) 



(31) 



(32) 



1 



Here we identified 

9{r) 



exp 



{x{t)x{t + t)) ^ 



-ijdr. 
(33) 

(34) 



Expressions {x{t)x{t + T))st and {x^)st are calculated in 
App.O 



C. Persistent and antipersistent cases 

Consider now the asymptotic behavior of the integrand 
in expression p3p at r — > oo. 



{...} 



[1+5- 

(a;2 



-9{r). 



9{r) 



{^^)st^+9{r) 



- 1 



(35) 



I St 

Since ^(r) ^ r^^^^, the integrand decays slowly; the 
integral in Eq. itself converges for H < 1/2 and 
diverges for H > 1/2. 

Focusing at first on the antipersistent case we notice 
that according to Eq. (1551) the main contribution comes 
from the integrand estimated at g^r) ~ 1, which imme- 
diately leads to 



T ~ exp 



1 



(36) 



being a kind of generalization of the standard transition- 
state arguments to the FBM case. Recalling that for our 
harmonic potential, {x^) = DT[2H + 1), we obtain an 
estimate for the coefficient b{H) in the empirical formula 
for the escape time, Eq. ([T^ . Namely, we find 



b{H) = 



1 



T{2H + 1) 



(37) 



Eq. p7p provides a surprisingly good approximation 
to the behavior of h{H) obtained from the simula- 
tions, as shown in Fig. [S] In particular, approximation 
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Figure 8: Coefficient b(H) occurring in the empirical formula 
(|12|) for the mean escape time. Symbols: Values obtained 
from best fit. Solid line: Theoretical behavior described by 
Eq. (|37|. 

P7p shows the nontrivial maximum for intermediate H- 
values. Fig. ini shows the values for the mean escape time 
obtained from our simulations of the antipersistent pro- 
cess with < H < 1/2, along with the behavior predicted 
by Eqs. §^ and 

In the persistent case the integral in expression ([33|) 
diverges. We show that a suitable truncation at some 
upper bound Tcut leads to a quite good agreement with 
the behavior recovered from simulations. Physically such 
a truncation always exists due to the finiteness of the 
slow power-law decay of the autocorrelation function of 
fractional Gaussian noise. Thus, we would always expect 
finite mean escape times also in the persistent range. Be- 
cause of the slow divergence of the integral for T, we 
may expect a weak dependence of the integral on the 
cutoff parameter Tcut if only it is chosen large enough. 
Indeed, we found in our numerical simulations that the 
value Tcut = 18 already gives good agreement with the 
numerical simulation, see Fig. IH] bottom. 



IV. SUMMARY 

In this work we present an extensive analysis of the 
generalized Kramers escape from a potential well for a 
particle subject to fractional Brownian motion. Specific- 
ally we considered a particle whose motion is governed 
by the Langevin equation driven by external fractional 
Gaussian noise. The motion we consider is thus not sub- 
ject to the fluctuation dissipation theorem. Potential ap- 
plications for such behavior may, for instance, include 
geo- and astrophysical fluctuations, stock market pricing, 
or teletrafRc. 

Based on simulations and analytical derivations we 
showed that, despite the driving fractional Gaussian 
noise, the escape dynamics preserved the classical expo- 
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Figure 9: (Color online.) Mean escape time T as function 
of inverse noise intensity, 1/D. Symbols: simulations results. 
The solid lines show the analytical result given by Eqs. 
and p4p in the antipersistent case < H < 1/2. For per- 
sistent motion 1/2 < ii" < 1 the solid lines represent a fit 
by Eqs. ((33)l and (|34|) based on numerical truncation of the 
integral in Eq. (|33|l with cutoff time Tcut = 18. 



nential shape of the distribution of escape times. Devi- 
ations from the behavior for regular Gaussian white noise 
are found in the activation dependence of the mean es- 
cape time on the noise intensity at different values of the 
Hurst exponent H. 

The escape turns out to slow down for increasing value 
of the Hurst exponent. Thus in the persistent case 
1/2 < H < 1 the escape is slower than in the antipersist- 
ent case < H < 1/2, and the latter is faster than for 
ordinary Brownian case. This somewhat surprising result 
is in accordance with previous results for the first passage 
time Q , where the scaling exponent of the first passage 
time distribution decreases for increasing H. We note 
that this observation is not restricted to the asymmetric- 
ally truncated harmonic potential used in this work, but 
also occurs for a symmetric truncation of the harmonic 
potential at x — ±\/2. 

Analyzing the detailed behavior of the mean escape 
time we find that the logarithm, logT in the entire sim- 
ulations range H — 0.1, ... , 0.85 depends linearly on the 
inverse noise intensity, 1/D. This activation dependence 
is thus preserved for both antipersistent and persistent 
cases. Conversely, the behavior of logT on the Hurst ex- 
ponent shows a linear dependence in the persistent case, 
while in the antipersistent case we find a nonlinear de- 
pendence. 

We note that fractional Brownian motion is an ergodic 
process in the sense that time and ensemble averages co- 
incide, albeit the convergence to ergodicity is algebraic- 
ally slow with the measurement time [56j . For sufficiently 
long averaging times the dynamic behavior of time and 
ensemble averages of individual trajectories should there- 
fore be identical. This contrasts the behavior for continu- 
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ous time random walk processes with diverging charac- 
teristic waiting times ^57] or with correlations in waiting 



si. 



times or jump lengths 

The understanding of fractional Brownian motion in 
several aspects remains formidable. We expect that this 
work contributes toward the dcmystification of this seem- 
ingly simple stochastic process. 
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Appendix A: Description of FBM generators 

Here we briefly describe the generators with which we 
simulated FBM. It should be noted that the generators 
provide best results for either the antipersistent case < 

< 1/2 or for the persistent case 1/2 < H < 1. 

A fast and precise (see the tests in Appendix B) gen- 
erator for fractional Gaussian noise in the anti-persistent 
case is described in Ref. [ij. In brief, the idea is as 
follows. 

First, we define a function 



[1 - (n/7V„,a: 

i?^(2Afmax-n), 



, , , , .x)'^J , < n < 7V„,ax 

^max <n< 27Vi„ax 

(Al) 

where H is the Hurst parameter (0 < < 1/2), n is the 
number of steps corresponding to time in the continuous 
time limit, and A^max is the length of the random sample. 
Second, we perform a discrete Fourier transformation of 
Eq. ini, with S^ik) = F{R^{n)}. 
We then define 



0, 

exp(i0fe)e(fc)v/^(fcy. 



X*(27V,„ax-fc), 



fc = 

< fc < iV,„ax 

k — N 

1^ — ^ max 



where the symbol * stands for complex conjugation. Ok 
are uniform random numbers from [0,27r), and ^(fc) are 
Gaussian random variables with zero mean and variance 
equal to 2. All random variables are independent of each 
other. 

Finally, we set y{n) — x{n) — a::(0), where x{n) — 
{X{k)} is the inverse Fourier transformation of 
Eq. (IA2I) . The quantity y{n) represents a free [i.e., in 
absence of an external force] fractional Brownian traject- 
ory which is to be differentiated with respect to time, to 



obtain fractional Gaussian random numbers. Since the 
variance depends on the number of steps A'max, it is 
normalized such that (^^) — 2. 

Despite the availability of several exact simulation 
methods, for the persistent case we chose an approximate 
but eflficient simulation method. This generator exploits 
the spectral properties of fractional Gaussian noise ji^ . 
The method uses the following steps: 

(i) Take white Gaussian noise ^(t), where t is an in- 
teger. 

(ii) Calculate the spectral density of this Gaussian 
noise and perform a Fourier transformation, S{k) — 

nm- 

(iii) Introduce correlations multiplying it by 1/k^ ^Z^, 
where 1/2 < iJ < 1. 

(iv) Inverse Fourier transform ^//(i) = 
F-i{5(fc)fci/^-^}, to obtain approximate fractional 
Gaussian noise with the index H. 

(v) Normalize the noise. 

In Appendix B we demonstrate that this method reli- 
ably produces FBM. 

We note that since we approximate the integral rep- 
resentation, this creates two types of errors, a 'low fre- 
quency' one due to the truncation of the limit of integra- 
tion and a 'high frequency' one caused by replacing the 
integral by a sum. By using various tests, we estimated 
the best discretization parameters. We used the max- 
imum sample length of 2^'* w 1.7 x 10^ steps, the time 
increment varying within the interval [0.001,0.01]. 



Appendix B: Testing the numerical algorithm 

To check our simulations algorithm based on numerical 
integration of the Langevin equation ()lip we performed 
a number of tests to validate the FBM we create with 
the generators sketched in Appendix A. 

First, we calculated the autocorrelation function of the 
fractional Gaussian noise. As shown in Fig. 1101 the sim- 
ulated data show excellent agreement with the analyical 
result (solid lines) given by Eq. ([5]) for discrete time steps. 

Second, we calculated the position mean squared dis- 
placement 



{xHity 



2Dt 



2H 



(Bl) 



and two-point correlation function 

{xH{tl)xH{t2)) - D{tl" + tl" - |ti - t2p^). (B2) 

of free FBM, and compare with the analytical expressions 
for FBM in discrete time n with time increments 5t = 1, 

{xH{nf) = 2Dn^", (B3) 
{xH{n)xHil)) = I?(l + n2^-|n-lp^). (B4) 

As demonstrated in Figs. [TT] and [T^l respectively, the 
agreement is excellent. 
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n 



Figure 10: Absolute value of the autocorrelation function of 
fractional Gaussian noise for the entire range of the Hurst 
exponent H in a, log-log scale as function of the number of 
time steps n. Inset: autocorrelation functions with the same 
Hurst indices for small numbers of steps n on a linear scale. 
The numerical results are shown for H = 0.10, 0.45, 0.60, 
and 0.80. The solid lines in the main graph correspond to 
the analytical solution (|8]). Simulations parameters: number 
of simulated samples A'stat ~ 20, 000, each of length A'^max = 
2^^ 8.2 X 10^ for the antipersistent case, and 2^^" 3.3 x lO" 
for the persistent case, respectively. 




10" 10^ 10^ 10^ 



n 

Figure 11: Mean squared displacement for free FBM in log- log 
representation. The solid lines show the analytical expression 
(|B3|) while the symbols depict the simulations for different 
Hurst parameters ranging from _ff = 0.1 (lowest curve) to 
H = 0.8 (uppermost curve). Here, D was taken to be equal 
to 1, the time step St = 1, iVstat = 20, 000, and iVmax = 2^" ^ 
10^ 



Third, solving Eq. (|lll) wc calculated the mean squared 
displacement for a particle in an infinite harmonic poten- 
tial well, as shown in Fig. 1131 The initial condition was 
a; = 0, at the bottom of the potential well. The asymp- 
totic analytical behaviors are represented by the initial 
free behavior {x'^{t)) ~ t^^ and the terminal saturation 



Figure 12: Position autocorelation function of free fBm. The 
solid lines show the analytical expression (IB4|) while the sym- 
bols depict the simulations for different Hurst parameters ran- 
ging from H = 0.1 (lowest curve) to H — 0.8 (uppermost 
curve). Again, D was taken to be equal to 1, the time incre- 
ment St^l, Astat = 2 X 10^ and A^ax = 64. 

value = Dr{l + 2H) at t ^ (X) (for details, 

see Appendix [n]) . This demonstrates that our generators 
also produce reliable behavior in an external potential. 

Finally, we performed a simulation of a free particle 
escaping from a semi-infinite axis with absorbing bound- 
ary under the influence of fractional Gaussian noise, see 
Fig. [TH The observed scaling of the first passage time 
PDF p (t) compares very favourably with the analytical 
result from Refs. [33 |: 

p{t) ~ ^25) 

Note that this relation cannot be obtained by the method 
of images, despite the fact that FBM is a Gaussian pro- 
cess. Also note that the slope of this power-law be- 
comes flatter for increasing Hurst coefficient: the escape 
is slower for a more persistent FBM, i.e., a motion whose 
mean squared displacement grows faster. This a priori 
surprising behavior is also seen for the escape from the 
potential well studied herein. 



Appendix C: Variance and autocorrelation function 
for FBM in a harmonic potential well. 

We now consider FBM in a harmonic potential, as 
described by the Langevin equation (compare with 

Eq. (nni)) 

^^-ax + D'/'^H{t), (CI) 

where we introduce the prefactor a which allows us to 
consider the harmonic potential {a — 1) and a free FBM 
(a = 0) as well. The solution of Eq. (|C1I) with the initial 
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Figure 13: Mean squared displacement of FBM in an har- 
monic potential. Top panel: mean squared displacement for 
fixed D = 1.0 and varying Hurst index. Bottom panel: mean 
squared displacement for fixed Hurst index H = 0.25 and 
four different values of the noise intensity. The solid lines 
show the asymptotes of free FBM {t^^ power-laws corres- 
ponding to straight lines in the log-log scale) and stationary 
states (horizontal lines) . Insets: stationary values of the mean 
squared displacements as functions of H for fixed D = 1.0. 
The points in all graphs represent the simulations results for 
the following parameters: time increment St = 0.01, num- 
ber of samples A'^stat = 10®, and number of steps per sample 
A'^max ~ 2^^' ~ 10'' for both persistent and antipersistent cases. 



condition a::(t = 0) = is 

ft 

x{t) = D^'^ / e-°(*-*')Cff(i')'^<'- (C2) 
Then, the ACF function 

tl t2 

{x{ti)x{t2)) ^ De-^'^' j dt' J dt"e^^''+'"'> iat'Mt")) 



tl t2 2 

= _7^e-2-* J dt' J dt"e''^''+*"^ \t' - i'f ?C3) 






o 
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Figure 14: First passage time PDF of FBM on a semi- 
infinite axis with absorbing boundary condition. The solid 
lines demonstrate the respective analytical slopes. Paramet- 
ers in the antipersistent case: H = 0.25, time increment 
5t — 0.001, number of samples A'stat = 100,000 and num- 
ber of steps per sample A'max ~ 2^^ ~ 1.3 x 10^. In the 
persistent case we used H = 0.75, St = 0.001, iVstat = 20, 000, 
and A^max = 2" « 8.2 x 10^ 



Now, if t2 — ti = T > 0, after some lengthy calculations 
we get Eq. jnH: 
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{x{t)x{t + T)) = i?|e-"('+^)i2« + e-'^*(t + T)2«-r2» 
2a2- 1 



2a(2i7 + 1) L^'''+'e-«(2*+r)^^(2if + l;2H + 2; at) 
+ {t + Tf"+^e-''^'^'+^'>M{2H +1;2H + 2; a{t + r)) 

- T^"+^e~''''M{2H + 1:2H + 2;aT) 

- ia-2(^+i) (2a2 _ 1) U^(r(2iJ + 1; a(f + t)) - r(2H + 1; ar)) 



e-''''(T(2H + l;at) - T{2H + 1)) 



(C4) 



Assuming a — 1, 



{x{t)x{t + T)) = Dle-'-'+^h^" -T^" + e-*{t + T) 



2H 



1 
1 



T{2H + 1) - T{2H + 1; <) + — — -M{2H + 1; 2i7 + 2; t) 

zii + 1 



(C5) 



'2t-T 



Af (2i/ + 1; 2i/ + 2; t) + ^ ' M{2H + l;2H + 2;t + t) 



2H + 1 



2H+1 



-e^ [T{2H + 1; r) - T{2H + 1; t + r) 



r 



Here, T{a, b) is the incomplete F-function, and M de- The autocorrelation function of free fBm can be nat 

notes the Kummer function [s^. In the stationary state urally obtained by placing a = in Eq. (|C4p : 
(t oo) the autocorrelation function Eq. (jCSP yields 



\2h 



{xit)x{t + T)) ^ D {t^'' + {t + r) 
{x{t)x{t + T))^f = dI [e-^r{2H + 1) + e^r(2i? + 1, r)] that matches the well-known relation M 



2h 



} 



(C8) 



2{2H + 1] 



M{2H + 1;2H + 2,t)-t 



2H 



In order to obtain the variance we take r = in 
Eq. (103: 

(x^(t)) = i:)i2t2^e"* + [r(i + 2i/)-r(i + 2i/,t)] 



(C6) 

Appendix D: Mean escape time and first escape 
time PDF for harmonic potential truncated from 
both sides 



2H+1 



2H + 1 



e-^^ M{2H + l;2H + 2]t) 



Now, the stationary variance is: 

{x^)^^ = DV{2H + l). 



(C7) 



Note that at r = Eq. (|C6l) reduces to Eq. ((C7l) . 
whereas for H = 1/2 it gives {x{t)x{t + t)) ^^ = De~'^ , 
the autocorrelation function of the Ornstein-Uhlenbeck 
process. Taking the asymptotics of the incomplete F- 
function and the Kummer function, one may easily see 

^2H-2 



In this Appendix we consider the Kramers problem 
for an harmonic potential, but this time we introduce a 
cutoff on both sides, that is, at a; = ±\/2, and evaluate 
the same dependencies (see Figures [T5l and [TO]) . 

One can observe that qualitatively there is no differ- 
ence in behaviour with the case of the one-side truncated 
potential. Indeed, the escape is faster when lowering the 
Hurst parameter; the escape time PDF remains exponen- 
tial and so does the mean escape time. Again, the MET 
may be fitted with the following function: 



T{H < 0.5) 
T{H > 0.5) 



exp(aa;^ + bx + c) 

exp(5'a; + c'), (Dl) 



that {x{t)x{t + t)) « 2DH{2H — 1)t at t — > oo. where a, 5, c, 6', c' are some constants depending on D. 
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Figure 15: First escape time PDF for harmonic potential trun- 
cated from both sides. Points are the simulation data, solid 
lines stand for linear fitting. Simulation details are the fol- 
lowing: for the antipersistent case D = 0.25, 5t — 0.001, 
A^max ~ 131072, A^stat = 10^; for the persistent case D = 
0.25, St = 0.002, Af,^ax = 131072, A^.tat = lO". 
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Figure 16: Mean escape time as function of the Hurst para- 
meter for harmonic potential truncated from both sides. 
Points are the simulation data, solid lines stand for fitting 
with Eq. (IDlf) . Simulation details are the following: for both 
antipersistent and persistent cases St varied from 0.001 to 
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